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Illustrations 


Figure  1.  Electron  (dot-dashed  line)  and  ion  densities  dotted,  He^ 
dashed,  O'*"  too  small  to  be  seen)  at  the  apex  of  the  L=3  fieldline  for  equinox, 
solax  moderate  conditions  plotted  as  functions  of  the  logarithm  of  density 
of  artificially  added  neutral  species,  demonstrating  reasonable  insensitivity 
to  that  parameter.  Panels  side-by-side  show  the  same  results  for  different 
methods  of  linearizing  the  equation  of  continuity.  Rows  of  panels  correspond 
to  results  for  0200LT,  lOOOLT  and  2200LT  respectively,  read  from  the  top. 

Figure  2.  As  in  Figure  1,  but  for  densities  at  the  apex  of  L=3.8. 

Figure  3.  A  comparison  of  modeled  GTIM  0"^  densities  with 
observed  profiles  from  the  Arecibo,  Puerto  Rico,  ISR.  The  four  panels  show 
comparisons  at  1900LT,  2200LT,  OlOOLT  and  0400LT.  Observations  are 
indicated  by  25  April  1995  as  solid  lines,  26  April  as  dotted  lines  and  27  April 
as  dashed  lines.  Overplots  show  radar  observations  made  within  15  minutes. 
Symbols  denote  model  densities  for  a  4-day  simulation  (day  l=plus,  day 
2=asterisk,  day  3=diamond,  day  4=triangle). 

Figure  4.  As  in  Figure  3,  but  pertaining  to  densities.  The 
additional  symbols,  “G”  demonstrate  the  results  of  the  Gallagher  empirical 
plasmaspheric  model. 

Figure  5.  As  in  Figure  3,  but  pertaining  to  He^  densities. 

Figure  6.  The  two  panels  refer  to  the  seasonal  variations  of  Nmax  at 
1900  and  2100LT  respectively.  Each  panel  shows  plots  of  Nmax  at  selected 
latitudes  as  a  function  of  month  during  the  year  at  the  US  sector.  While  the 
results  at  ±10°  (dotted  and  dashed  lines)  are  typically  very  similar,  those  at 
±20°  have  been  as  a  measure  of  the  asymmetry,  with  the  hand-drawn  insert 
used  to  reflect  the  level  of  symmetry. 

Figure  7.  As  in  Figure  6,  but  for  the  Pacific  sector. 


Figure  8.  GTIM  results  for  Nmax  and  Hmax  plotted  against  Millstone  20 

Hill  radar  measurements.  The  model  has  been  run  with  both  climatological 
neutral  winds  and  with  measured  winds.  The  model  run  with  adjusted  winds 
show  significant  improvement  in  Hmax  at  night.  The  enhanced  equatorward 
winds  were  able  to  maintain  the  nighttime  ionosphere  (Nmax)  under  active 
magnetic  conditions. 

Figure  9.  Time  variations  of  modeled  Nmax  values  before  and  after  the  23 

onset  of  a  magnetic  disturbance,  at  five  selected  latitudes.  Leftmost  in  each 
plot  is  the  quiet-time  variations,  and  this  is  repeated  as  a  dashed  line  for  the 
subsequent  days.  The  solid  line  shows  the  variation  in  Nmax  following  storm 
onset,  and  positive  and  negative  storm  phases  being  clearly  apparent. 

Figure  10.  Wavelike  ionospheric  effects  in  the  GTIM  when  a  wavelike  E  26 

X  B  disturbance  is  applied.  This  case  was  made  for  January,  solar  maximum, 

Huancayo  sector  conditions,  using  a  16-day  wave  and  50%  peak-to-peak 
amplitude.  The  relative  variation  in  Nmax  at  selected  latitudes  (18=solid 
line,-10=dotted,0=dashed,-fl0=dot-dashed,-|-18=dot-dot-dot-dashed)  is 
plotted  as  a  function  of  the  day  during  the  period. 

Figure  11.  A  demonstration  of  the  abilities  of  the  SIMP  LI  profile  model  29 

in  radiowave  propagation  predictions.  Here,  the  solid  line  indicates  the 
elevation  angle  versus  transmission  frequency  for  a  given  circuit  obtained 
from  climatological  estimates  (electron  density  contours  in  the  left  insert). 

The  dotted  line  shows  the  same  predictions  made  from  the  “ground-truth” 
ionosphere,  shown  in  the  right  insert.  The  three  topmost  lines  show  the 
results  when  Nmax  is  updated  at  the  circuit  mid-point,  end-points,  and  both, 
and  when  SIMPLI  parameters  along  the  path  are  shifted  only  the  zeroeth 
order.  The  three  lines  lying  nearest  to  the  dotted  line  correspond  to  Nmax 
being  refit  along  the  whole  path. 

Figure  12.  Averaged  periodogram  for  quiet  vs.  active  orbits  in  1983.  37 
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Figure  13.  Day  vs.  night  averaged  periodogram  for  1983.  38 

Figure  14.  Kp  vaxiations  for  the  July  ,  1983,  and  August,  1983  study  42 

periods. 

Figure  15.  Daytime  density  variations  during  July  20-27,  1983.  43 

Figure  16.  Same  as  Figure  15,  except  for  nighttime.  44 

Figure  17.  Same  as  Figure  15,  except  for  August,  1983  study  period.  45 

Figure  18.  Comparison  of  SETA  and  MSISE90  density  variations  47 

averaged  over  20-20  degrees  latitude  at  1030  LT  during  July  20-26, 1983. 

Figure  19.  Comparison  of  ”longitude/UT”  variation  in  SETA  densities  49 


(circles)  with  MSISE90  (dotted  line)  during  the  magnetically  quiet  period  of 
July  21-23,  1983.  Solid  line  represents  least-squares  fit  of  sinusoid  to  SETA 
data. 

Figure  20.  Comparison  of  ampltiudes  of  ”longitude/UT”  variation  50 

between  SETA  and  MSISE90  as  a  function  of  latitude,  for  daytime  data 
during  July  21-23,  1983. 

Figure  21.  Longitude/UT  dependence  as  revealed  by  SETA  (top)  and  51 

MSISE90  (bottom)  total  mass  densities  at  200  km  during  July  21-23,  1983. 

Left:  daytime  («1030  LT).  Right:  nighttime  («2230  LT). 

Figure  22.  Longitude/UT  dependence  as  revealed  in  SETA  (left)  and  52 

MSISE90  (right)  daytime  (wl030  LT)  total  mass  densities  at  200  km  during 
December  5-10, 1983. 

Figure  23.  Density  residuals  as  a  function  of  latitude  and  time  during  55 

magnetically  quiet  periods  during  July,  1983,  and  December,  1983. 

Figure  24.  Superposition  of  spectra  in  different  latitude  bands,  56 

corresponding  to  the  density  residuals  depicted  in  Figure  17  for  July,  1983. 
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1  Introduction 


In  order  to  provide  a  better  description  of  the  ionospheric  environnaent  a  variety  of 
studies  have  been  undertaken  covering  a  range  of  different  approaches.  These  span  the 
spectrum  from  the  expansion  of  scope  of  a  theoretical  ionospheric  model,  passing  through 
a  series  of  event-oriented  modeling  studies,  to  the  investigation  of  new  applications-oriented 
empirical  descriptions  of  ionospheric  behavior. 

From  the  empirically-based  applications  studies,  new  quick  methods  can  be 
provided  to  users  of  real-time  and  monitoring  systems  to  enhance  their  performance 
and  capabilities,  and  new  clues  to  underlying  physical  processes  may  be  gleaned.  From 
the  event-oriented  modeling  studies,  the  validity  of  physical  interpretation  of  particular 
ionospheric  phenomena  can  be  determined.  From  the  expansion  of  a  theoretical  model  to 
multiple  ions  and  much  greater  altitudes,  the  scope  of  investigations  are  greatly  enhanced 
and  a  new  generation  of  specification  model  can  ultimately  be  provided  to  the  user. 

The  Phillips  Laboratory  Global  Theoretical  Ionospheric  Model  (hereafter,  GTIM) 
has  been  rewritten  to  solve  for  multiple  ion  densities,  and  physical  drivers  have  been 
updated.  The  GTIM  has  been  utilized  to  demonstrate  (1)  that  latitudinal  symmetry 
in  equatorial  electron  densities  can  be  related  to  scintillation  occurrence  rates  (2)  that 
ionospheric  responses  to  substorms  can  be  investigated,  (3)  that  neutral  winds  can  play  an 
important  role  in  nighttime  maintenance  and  (4)  that  satellite  in  situ  measurements  can 
be  used  to  derive  F-region  parameters.  The  model  has  further  been  adapted  to  the  study 
of  disturbed  ionospheric  conditions,  both  of  impulsive  and  periodic  types.  Finally,  on  the 
applications  side,  studies  were  conducted  showing  (1)  that  a  new  simple  profile  model  has 
powerful  ray-tracing  applications,  (2)  that  the  same  model  can  be  used  to  combine  satellite 
measurements  from  different  sources  and  (3)  that  GPS  satellites  offer  an  opportunity  for 
monitoring  scintillations. 
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On  the  thermospheric  side,  Satellite  Electrostatic  Triaxial  Accelerometer  Experiment 
(hereafter  SETA)  density  measurements  have  been  extensively  utilized  to  both  investigate 
and  subsequently  describe  a  series  of  phenomena  (1)  meridional  wave  structures,  (2) 
magnetic  storm  response  and  (3)  longitudinal  wave  structures. 

2  Extending  the  GTIM  to  the  Plasmasphere 

One  significant  effort  over  this  period  has  been  the  extension  of  the  Phillips 
Laboratory  GTIM  to  the  solution  of  multiple  ion  species  over  a  significantly  greater  range 
of  altitudes.  This  has  involved  two  major  areas  of  development  within  the  model  that  were 
performed  at  Boston  University: 

-  rewriting  the  equations  of  continuity  and  momentum  for  more  than  one  ion, 
including  both  ion-ion  and  ion-electron  interactions 

-  updating  physical  drivers  for  the  model  to  permit  reliable  ion  density  derivations  at 
much  higher  altitudes. 

There  has  also  been  a  significant  effort  made  in  validating  the  model  results  and  in 
verifying  the  numerical  stability  of  the  ion  densities.  The  multiple-ion  GTIM,  as  described 
briefly  in  the  following  section,  is  considered  to  be  reliable  up  to  L=3. 

Working  code,  and  revised  model  drivers  have  been  delivered  to  colleagues  at 
Computational  Physics  Inc.,  in  Boston.  The  multiple  ion  model  was  described  in  a 
presentation  at  the  1996  AGU  Meeting,  in  San  Francisco,  CA,  and  is  currently  in 
preparation  as  a  journal  article. 

2.1  Deriving  Ion  Densities 

The  GTIM  solution  method  for  a  one-ion  formulation  is  described  in  Anderson  (1973). 
Essentially,  an  expression  for  ion  flux  derived  from  the  momentum  equation  is  substituted 
into  the  continuity  equation.  This  expanded  continuity  equation  is  then  finite-differenced. 


2 


and  the  resulting  equation  that  is  tridiagonal  in  the  ion  densities  (denoted  Ni),  is  solved 
numerically  along  magnetic  field  lines  from  the  footpoint  at  125km  in  one  hemisphere  to  the 
conjugate  point  in  the  other  hemisphere. 

We  now  consider  the  multiple  ion  situation.  In  the  first  instance,  we  form  an 
equivalent  expression  for  ion  fluxes  from  the  momentum  equation  for  each  of  a  number  of 
ions.  The  momentum  equations  include  terms  from  each  of  collisions  with  other  ions  (say,  j) 
and  with  electrons.  The  electron  momentum  equation  is  used  to  eliminate  the  electric  field 
term  in  the  manner  of  the  one-ion  formulation,  and  the  resultant  expression  for  ion  fluxes  is 
of  the  same  general  form  as  that  obtained  for  one  ion,  permitting  the  same  substitution  into 
the  continuity  equation  and  the  same  type  of  numerical  solution.  However,  in  three  ways, 
the  expression  is  more  complex  than  the  one-ion  case: 

-  ion  density  and  electron  density  are  no  longer  equal,  and  NijNe  terms  appear  in  the 
expression 

-  ion-ion  interactions  axe  included,  and  ion  velocities  become  important  to  the  solution 
as  the  relative  velocities  between  ions  define  the  momentum  transfer 

-  ion-ion  couplings  for  the  ith  ion  also  include  coupling  of  the  jth  ion  with  electrons, 
owing  to  the  electric  field  term  in  the  momentum  equation. 

Thus,  strictly  speaking,  the  equation  for  ion  flux  is  not  linear  in  Ni.  On  the  other  hand,  the 
higher  order  terms  are  generally  small  and  can  be  arranged  to  appear  only  as  modifiers  to 
collision  frequencies.  Rearranging  the  higher  order  terms  and  determining  these  effective 
collision  frequencies  are  key  steps  in  obtaining  density  solutions  in  the  multiple  ion  GTIM. 

We  have  thus  outlined  an  approach  for  deriving  ion  densities  for  each  of  a  number  of 
different  species.  The  ions  required  to  define  the  ionosphere  and  plasmasphere  are  0"^, 
and  and  solutions  are  sought  for  each  ion  in  turn  at  each  time  step.  Because  the  ion 
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densities  are  solved  sequentially  rather  than  simultaneously,  much  care  is  taken  to  ensure 
that  ions  are  not  artificially  created  or  destroyed  as  numerical  artifacts  in  any  time  step. 
This  involves  making  adjustments  to  both  production  and  loss  rates  and  is  best  illustrated 
in  an  example. 

Let  us  look  at  the  charge  exchange  reaction  and  first  consider  the  situation  at 

an  altitude  of  around  300km.  The  loss  rate  (/3)  for  O'*"  due  to  this  charge  exchange  amounts 
to  typically  30%  of  the  total  O'*"  loss  rate,  and  thus,  when  this  reaction  is  included  in  the 
calculations,  model  F-region  O'*"  densities  drop  by  around  30%.  However,  at  these  low 
altitudes,  no  if"*"  can  be  formed,  owing  to  the  rapid  loss  of  if"*"  through  the  reverse  reaction. 
In  this  situation,  the  answer  is  to  modify  the  forward  loss  rate  by  a  factor  that  depends  on 
the  reverse  loss  rate  {^rev)  and  the  time  step  (AT).  We  say  that  the  fraction  of  O'*"  that 
can  be  lost  out  of  the  amount  one  is  trying  to  lose  is  the  same  as  the  fraction  of  that  one 
is  able  to  form  over  the  time  step,  and  that  the  rest  charge  exchanges  rapidly  back  to  O"*". 
Specifically,  this  fraction  is  exp{-^REV^T). 

At  a  higher  altitude  (say,  1500km),  is  primarily  a  source  of  O'*".  One  derives  a 

production  term  (locally)  for  0"*"  that  goes  as  k  n(jH'*')  n(0)  AT.  However,  this  production 
should  not  physically  exceed  the  density  of  the  available  source  terms  (here  jH"*"  and  0, 
mostly  the  latter  being  the  smaller  term  at  high  altitudes).  To  ensure  this,  production 
rates  are  modified  by  a  factor  of  exp(-«  AT)  to  reflect  the  fraction  of  0  that  can 

be  lost  over  the  time  step.  It  should  also  be  noted  that  because  both  n(0)  and  ^REviH"^) 
are  steeply  varying  functions  of  altitude,  that  the  altitude  profiles  of  these  exp(-X  AT) 
terms  are  effectively  step  functions,  and  that  the  choice  of  how  one  averages  these  decay 
terms  over  the  time  step  (exp(-X  AT)  vs  exp(-X  AT/2),  say)  has  no  influence  on  the  final 
densities  derived  in  the  model. 
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2.2  Revising  Model  Drivers 


As  outlined  above,  each  ion  density  can  be  derived  from  an  expanded  continuity 
equation  (using  effective  ion  fluxes),  that  becomes  a  tridiagonal  matrix  when  finite- 
differenced.  The  coefficients  that  relate  the  adjacent  terms  in  this  matrix  are  dependent  on 
a  number  of  factors.  Some  of  these  factors  are  those  that  define  the  physical  environment 
in  which  the  ion  density  solutions  are  being  derived,  such  as  thermospheric  densities  and 
winds,  electron  and  ion  temperatures,  and  vertical  E  x  B  drifts.  In  the  GTIM  formulation, 
these  physical  quantities  are  not  derived  but  are  simply  input,  and  are  referred  to  as 
physical  drivers  for  the  (ionospheric)  model. 

Two  important  drivers  have  undergone  recent  revisions,  specifically  to  address  the 
issue  of  the  higher  altitudes  being  studied  in  the  multiple  ion  GTIM. 

First,  vertical  drifts  have  been  generalized  to  include  a  vaxiation  with  altitude. 

Most  commonly  this  takes  the  form  of  using  input  values  (based  on  Jicamarca  or  AE 
observations)  at  low  altitudes,  tapering  off  above  one  input  altitude  to  zero  by  a  second 
input  altitude.  An  additional  term  in  the  continuity  equation  resulting  from  the  non-zero 
dV/dZ  has  been  included.  Furthermore,  because  this  term  can  influence  ion  densities  at  the 
10%  level,  the  altitude  dependence  has  been  made  smooth  in  the  first  derivative  to  avoid 
discontinuities  in  densities  with  altitude  (runs  at  adjacent  L-shells).  A  tanh  dependence  is 
used.  The  input  altitudes  that  define  the  vertical  dependence  were  seen  to  vary  with  solar 
activity,  and  possibly  other  factors. 

Next,  the  specification  of  electron  and  ion  temperatures  has  been  revised.  Brace  and 
Theis  (1981)  provided  a  convenient  tabulation  of  a  summary  of  measured  Te  at  several 
altitudes  (300km,  400km,  1400km  and  3000km).  At  each  location  in  the  GTIM,  a  system 
has  been  set  up  that  provides  a  smooth  fit  with  altitude  to  the  Brace  and  Theis  values  at 
that  location,  using  a  tanh  functional  form.  The  fit  is  further  constrained  by  the  following 
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conditions: 


-  Te  is  set  to  the  neutral  temperature  (from  MSIS)  at  the  fieldline  baise  altitude 
(125km) 

-  Te  being  monotonic  increasing  with  altitude  at  any  location  (not  the  same  as 
monotonic  increasing  along  the  fieldline) 

-  Te  being  smoothly  and  slowly  varying  (larger  rates  of  change  lead  to  anomalously 
high  ion  velocities  in  the  model  that  in  turn  lead  to  large  density  gradients  and  sometimes 
numerical  instabilities). 

The  second  and  third  points  mean  that  the  fit  will  not  necessarily  pass  through  the  Brace 
and  Theis  values,  but  will  instead  be  a  smooth  fit  guided  by  those  values. 

At  higher  altitudes  than  3000km,  it  is  not  sufficient  to  simply  set  Te  to  the  3000km 
value.  In  filling  the  plasmasphere,  the  model  is  basically  diffusing  light  ions  up  both 
hemispheres  of  a  fluxtube  until  the  pressure  gradient  term  balances  the  gravity  term  (in 
the  momentum  equation).  A  flat  Te  profile  above  3000km  would  mean  that  much  larger 
light  ion  densities  are  required  at  high  altitudes  to  achieve  this  balance,  and  that  the  model 
refilling  time  is  very  long  (tens  of  days).  Above  3000km,  we  have  included  a  simple  linear 
gradient  in  Te  with  altitude,  and  to  avoid  ion  velocity  related  instabilities  in  the  solution 
from  steep  temperature  gradients,  there  is  an  additional  transition  region  above  3000km 
where  the  Te  gradient  changes  linearly  with  altitude  from  the  3000km  value  to  the  input 
asymptotic  value  at  a  higher  altitude  (say,  4000km).  In  this  fashion,  light  ion  densities 
converge  within  ten  days  or  so,  depending  on  altitude. 

The  ion  temperature,  Ti,  is  derived  from  Te  and  via  a  simple  formulation  of  Banks 
and  Kockarts  (1973)  that  ensures  a  match  with  Tjy  at  low  ionization  fraction,  and  Te  at 
higher  altitudes/fractions. 
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2.3  Testing  Model  Stability  and  Reliability 


A  couple  of  important  issues  regarding  model  stability  and  reliability  have  been 
recently  addressed. 

The  first  issue  arises  from  the  higher  order  terms  in  Ni  that  appear  in  the  multiple  ion 
momentum  equation  for  the  ith  ion.  There  are  different  methods  of  rearranging  the  higher 
order  terms  to  make  the  expression  apparently  linear  in  Ni  to  facilitate  the  GTIM  solution 
for  densities  in  the  expanded  continuity  equation.  Thus,  it  is  important  to  show  that  the 
derived  densities  do  not  depend  on  the  method  of  rearranging  used. 

The  second  issue  arises  from  the  use  of  ion  fluxes  from  the  momentum  equation  in 
the  continuity  equation,  as  it  relies  on  collisions  between  species  to  define  those  fluxes.  At 
sufficiently  high  altitudes,  collisions  are  too  rare,  and  artificial  neutral  atoms  are  added 
to  achieve  numerical  stability.  Thus,  it  is  important  to  chronicle  when  these  additional 
neutrals  are  required,  and  when  they  start  to  influence  the  densities  derived. 

Figure  1  demonstrates  results  obtained  for  ion  densities  at  the  apex  of  the  L=3 
fieldline  for  equinox,  solar  moderate  conditions  for  a  five  day  run.  Each  panel  shows 
densities  of  each  ion  (and  electrons)  at  a  particular  local  time  as  a  function  of  the  logarithm 
of  density  of  neutrals  added.  Flat  lines  thus  indicate  the  desired  insensitivity  of  the  results 
to  this  numerical  artifact.  Adjacent  panels  show  the  same  results  using  different  methods  of 
treating  higher  order  terms  in  Ni,  and  the  overall  consistency  of  ion  densities  (to  10%  or  so 
in  this  case)  demonstrates  insensitivity  to  the  method.  Other  rows  show  the  same  behavior 
at  other  local  times.  The  results  are  seen  in  Figure  1  to  be  encouraging  at  L=3.  However, 
by  L=3.8,  as  seen  in  Figure  2,  the  sensitivity  to  both  neutrals  and  method  is  much  more 
substantial  (say,  50%),  demonstrating  that  a  different  solution  approach  will  be  required  at 
the  highest  altitudes  in  the  plasmasphere. 

At  lower  altitudes,  reasonable  agreement  with  observed  ion  densities  has  been 
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Oenstly  Density  Density 
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Figure  1.  Electron  (dot-dashed  line)  and  ion  densities  (Ff"*"  dotted,  Ht^ 
dashed,  O'*"  too  small  to  be  seen)  at  the  apex  of  the  L=3  fieldline  for  equinox, 
solar  moderate  conditions  plotted  as  functions  of  the  logarithm  of  density 
of  artificially  added  neutral  species,  demonstrating  reasonable  insensitivity 
to  that  parameter.  Panels  side-by-side  show  the  same  results  for  different 
methods  of  linearizing  the  equation  of  continuity.  Rows  of  panels  correspond 
to  results  for  0200LT,  lOOOLT  «ind  2200LT  respectively,  read  from  the  top. 


8 


Density  Density  Density 


mOGE 

L-  3.80  IVO  2LT 


L-  3.80  IVO  2LT 


L-  3.80  IVO  2LT 


FUDGE 
L-  3.80  IV012LT 


FUDGE 

L-  3.80  IV012LT 


FUDGE 

t-  3.80  IV012LT 


FUDGE 

L*  3.80  IV022LT 


FUDGE 

Im  3.80  IV022LT 


FUDGE 

L-  3.80  IV022LT 


Figure  2.  As  in  Figure  1,  but  for  densities  at  the  apex  of  L=3.8. 
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obtained.  Figures  3,  4  and  5  show  comparisons  of  modeled  profiles  of  0'^,  and  He^ 
respectively,  for  selected  local  times  and  Arecibo  ISR  measured  profiles  from  April 
1995.  In  each  panel,  the  different  linestyles  demonstrate  the  variations/uncertainties 
in  the  observations,  while  the  symbols  show  the  model  results  for  consecutive  run  days. 
Overlapping  symbols  demonstrate  convergence  in  the  model  densities,  and  qualitative  and 
quantitative  agreement  with  the  observations  is  quite  reasonable  to  2000km.  Comparisons 
have  also  been  made  with  climatological  DMSP  in  situ  measurements,  kindly  made 
available  by  Dr.  Fred  Rich  at  Phillips  Laboratory,  and  again  the  results  are  generally 
encouraging,  especially  for  the  major  ion. 

Overall,  He^  has  not  been  seen  in  the  model  to  ever  be  the  major  ion,  but  it  is  often 
evident  in  the  10-30%  of  total  electron  densities.  H'^-0'^  transition  heights  in  the  model 
are  seen  at  altitudes  as  low  as  500-600km,  in  line  with  many  observations  made  at  solar 
minimum. 

3  Modeling  Studies  Using  the  GTIM 

In  terms  of  being  able  to  perform  diagnostic  modeling  studies,  an  advantage  of  the 
GTIM  is  that  the  thermospheric  and  electrodynamic  environment  of  the  ionosphere  is 
input,  rather  than  being  derived  in  a  self-consistent  manner.  It  is  thus  straightforward  to 
investigate  the  ionospheric  effects  of  varying  a  variety  of  important  physical  drivers,  many 
of  which  will  not  be  known  in  an  operational  sense.  Modeling  studies  of  this  type  can  thus 
help  to  isolate  physical  processes  by  their  ionospheric  signatures.  Several  such  studies  are 
described  in  the  following  section. 

3.1  The  Equatorial  Anomaly-Irregularities 

Ionospheric  irregularities  are  a  typical  occurrence  at  equatorial  latitudes,  and  the 
morphology  of  their  behavior  has  been  described  observationally  both  in  terms  of  season 
and  longitude.  A  theoretical  modeling  study  of  these  features  in  the  equatorial  anomaly 
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cm^~3 
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Fimre  3.  A  comparison  of  modeled  GTIM  0"^  densities  with 
observed  profiles  from  the  Arecibo,  Puerto  Rico,  ISR.  The  four  panels  show 
comparisons  at  1900LT,  2200LT,  OlOOLT  and  0400LT.  Observations  are 
indicated  by  25  April  1995  as  solid  lines,  26  April  as  dotted  lines  and  27  April 
as  dashed  lines.  Overplots  show  radar  observations  made  within  15  minutes. 
Symbols  denote  model  densities  for  a  4-day  simulation  (day  l=plus,  day 
2=asterisk,  day  3=diamond,  day  4=triangle). 
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region  weis  performed  with  a  particular  view  to  isolating  disturbance  effects.  The  GTIM 
is  driven  by  three  inputs  that  vary  with  geomagnetic  activity.  These  are  (1)  the  neutral 
atmosphere  (from  MSIS)  (2)  the  neutral  wind  field  (from  the  NASA  model)  and  (3)  the 
vertical  E  x  B  plasma  drifts.  It  was  next  necessary  to  derive  storm-time  E  x  B  drift 
patterns.  Alterations  were  made  to  the  drift  patterns  in  three  ways,  that  seemed  to  best 
summarize  the  changes  in  the  drift  patterns  (1)  at  the  pre-reversal  enhancement  (PRE) 

(2)  towards  the  daytime  peak  (0900-1400LT)  (3)  towards  the  nighttime  minimum  (0000- 
0500LT). 

Model  runs  were  made  to  compare  with  a  recent  quantitative  summary  of  the 
morphology  of  the  occurrence  of  scintillation  at  low-latitude  stations  at  four  distinct 
longitude  regions,  characterized  by  the  presence  of  stations  and  the  magnetic  declination. 
The  results  were  examined  with  a  view  to  quantifying  the  level  of  latitudinal  symmetry 
in  densities  in  the  local  ionosphere.  Broadly  speaking,  the  months  at  each  location  where 
the  symmetry  in  Ne  is  highest  match  the  months  of  highest  spread-F  occurrence  (and  thus 
scintillations).  Here,  we  refer  to  Figures  6  and  7  to  show  the  seasonal  variations  of  density 
asymmetry  in  the  model  in  the  US  and  Pacific  sectors,  respectively.  The  two  panels  refer 
to  results  at  1900  and  2100LT  respectively.  Each  panel  shows  plots  of  Nmax  at  selected 
latitudes  as  a  function  of  month  during  the  year  at  the  named  sector.  While  the  results 
at  ±10®  (dotted  and  dashed  lines)  are  typically  very  similar,  those  at  ±20®  have  been  as  a 
measure  of  the  asymmetry,  with  the  hand-drawn  insert  used  to  reflect  the  level  of  symmetry. 
Referring  to  the  descriptions  of  Aarons  (1993),  the  US  sector  results  resemble  those  seen  in 
Huancayo,  and  those  in  the  Pacific  sector  resemble  Kwajalein,  strongly  peaked  in  equinox 
months. 

While  the  mean  thermospheric  wind  model  of  Hedin  et  al.  (1988)  is  apparently 
sufficient  to  model  the  climatology  of  spread-F,  as  demonstrated  above,  case  studies  should 
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USA  Sector  :  1985  19LT  .C 


Lots  ;  20  10  -10  -20 

USA  Sector  :  1  985  21LT  .C 


Figure  6.  The  two  panels  refer  to  the  seasonal  variations  of  Nmax  at 
1900  and  2100LT  respectively.  Each  panel  shows  plots  of  Nmax  at  selected 
latitudes  as  a  function  of  month  during  the  year  at  the  US  sector.  While  the 
results  at  ilO*  (dotted  and  dashed  lines)  are  typically  very  similar,  those  at 
i20  have  been  as  a  measure  of  the  asymmetry,  with  the  hand-drawn  insert 
used  to  reflect  the  level  of  s)Tnmetry. 


be  based  upon  measurements  of  neutral  winds  whenever  these  are  available,  given  the  role 
of  winds  in  determining  symmetry  at  low  latitudes. 

In  the  manner  outlined  above  the  same  modeling  study  has  been  conducted  including 
disturbance  effects  on  the  vertical  drifts.  Broadly  speaking,  the  seasonal  variations 
of  symmetry  in  modeled  electron  density  are  dilferent  during  disturbed  periods.  The 
prediction  from  these  model  runs  is  that  the  likelihood  of  scintillations  during  storms  at 
these  two  longitude  sectors  will  be  less  peaked  at  equinoxes,  with  occurrence  at  the  June 
solstice  becoming  almost  as  likely  as  at  equinoxes.  This  does  not,  however,  factor  in  the 
higher  occurrence  rates  of  storms  during  equinox  months. 

3.2  Substorms 

Here,  the  GTIM  was  used  to  simulate  total  electron  content  (TEC)  disturbance 
events  observed  at  middle  and  lower  latitude  sites  near  75® W  and  7°E  longitudes  during 
substorms.  Within  this  longitudinal  range,  the  daytime  TEC  disturbances  show  a  pattern 
in  which  a  positive  response  can  be  followed  by  a  negative  response  or  depletion,  or 
diurnal  double  maxima  (DDM)  patterns  if  substorm  activities  axe  repeated.  The  required 
dayside  electric  fields  in  both  longitudes  showed  eastward  and  westward  perturbations 
during  substorm  expansion  and  subsidence  phases,  respectively.  This  is  suggestive  of 
magnetospheric  electric  field  penetration  and  subsequent  over  shielding  effects.  The 
morphological  features  of  the  required  electric  field  perturbations  near  dawn  and  dusk  were 
examined  to  derive  the  local  time  characteristics  of  magnetospheric  influence.  Large  scale 
traveling  atmospheric  disturbances  (TADs)  were  also  investigated  in  the  model. 

The  ionospheric  positive  response  to  substorms  has  long  been  suggested  to  be  due  to 
TADs  (historically  travelling  ionospheric  disturbances,  TIDs,  referring  to  the  ionospheric 
effects)  in  middle  latitudes.  There  is  a  time  delay  of  the  ionospheric  peak  height  behavior  at 
middle  latitudes  relative  to  substorm  activity  shown  at  high  latitudes,  attributed  to  TAD 
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propagation.  This  raises  the  possibility  of  other  dynamical  sources  which  can  rapidly  affect 
the  global  ionosphere,  such  as  electrodynamic  effects.  In  addition,  the  TAD  mechanism  for 
producing  observed  TEC  disturbances  also  needs  to  address  the  issue  of  TAD’s  appearing 
from  the  opposite  hemisphere  as  substorms  subside,  as  well  as  their  long  duration  and  very 
rapid  propagation  speeds. 

We  have  identified  a  consistent  daytime  ionospheric  disturbance  pattern  at  middle 
and  lower  latitudes  during  substorms;  namely,  an  ionospheric  TEC  positive  response 
corresponding  to  substorm  expansion  periods,  followed  by  a  depletion  or  negative  response 
as  the  substorm  subsides.  This  work  in  this  section  was  published  as  Pi  et  al.  (1995). 

3.3  The  Nighttime  Ionosphere 

To  investigate  the  issue  of  nighttime  and  topside  maintenance  of  the  ionosphere, 
the  GTIM  was  run  at  middle  latitudes  and  compared  with  incoherent  scatter  radar 
measurements.  A  request  for  an  ionospheric  experiment  was  made  by  Xiaoqing  Pi  and 
David  N.  Anderson  to  the  Atmospheric  Research  Group,  Haystack  Observatory,  to 
measure  electron  density  profiles,  meridional  neutral  winds,  and  electron  fluxes  under 
conditions  of  solar  minimum,  winter,  and  particularly  at  night.  Such  an  experiment 
was  conducted  three  times  during  December,  1994  and  January,  1995,  of  about  24 
hours  duration  each,  using  the  Millstone  Hill  incoherent  scatter  radar  (ISR)  and 
the  coordinated  Fabry-Perot  interferometer  (FPI)  operated  by  the  Atmospheric 
Research  Group.  The  calibrated  radar  data  are  accessible  via  internet  at  a  web  site 
(http://hyperion.haystack.edu/homepage.html,  datafile=mil950130g.004). 

The  GTIM  was  run  under  both  climatological  and  real  conditions  appropriate  to 
January  and  April,  1995.  The  Kp  index  of  January  30-31,  1995  show  that  geomagnetic 
conditions  were  active  during  the  experiment.  The  model  with  standard  input  (MSIS  and 
the  NASA  wind  model)  gave  results  that  matched  the  daytime  ionospheric  electron  density 
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profiles  quite  well,  but  underestimated  the  nighttime  values.  In  the  topside,  the  0"^  model 
underestimated  the  measurements  above  600  km. 

Some  model  input  adjustments  have  been  made  to  match  the  nighttime  ionospheric 
measurements.  One  adjustment  was  to  replace  the  model  winds  with  the  measured 
meridional  neutral  winds.  The  wind  measurements  were  made  using  both  ISR  and  FPI 
techniques,  and  two  measurements  showed  much  larger  nighttime  equatorwaxd  winds 
(300  m/s  or  more)  than  the  model,  though  the  two  measurements  differ  to  some  extent. 
These  strong  nighttime  equatorward  winds  may  have  been  due  to  magnetic  activity.  The 
winds  derived  from  the  radar  measurements  have  been  taken  for  the  model  run  due  to 
its  availability  in  both  day  and  night.  The  model  run  with  the  wind  adjustment  gave 
good  agreement  in  nighttime  Nmax,  but  showed  higher  Nmax  values  than  the  daytime 
measurements.  The  model  results  with  adjusted  winds  also  showed  significant  improvement 
in  Hmax  at  night.  This  is  seen  in  Figure  8,  where  the  two  panels  compare  Nmax  and  Hmax 
model  results  with  the  radar  observations.  These  model  runs  indicated  that  the  enhanced 
equatorward  winds  were  able  to  maintain  the  nighttime  ionosphere  under  active  magnetic 
conditions  at  this  subauroral  site.  Maintenance  of  higher  densities  in  the  topside  may  well 
be  provided  by  ions  and  by  charge  exchange  with  0  acting  as  a  source  for  0^. 

I 

This  work  was  presented  at  the  1995  CEDAR  workshop  in  Boulder,  CO. 

1 

3.4  Satellite  Signatures 

The  goal  of  this  project  was  to  determine  what  may  be  deduced  about  the  dayside 
ionosphere  and  the  ambient  neutral  winds  and  vertical  drifts  given  only  an  in  situ  electron 
density  measurement,  Ne(840),  in  order  to  enhance  the  real-time  ionospheric  monitoring 
capabilities  of  the  DMSP  satellite  system. 

A  modeling  study  was  undertaken  using  the  O'*"  GTIM.  Important  inputs  to  the 
model  are  the  specified  diurnal  variations  of  both  the  vertical  ion  drifts  and  the  meridional 
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TT-n  GTIM  results  for  Nmax  and  Hmax  plotted  against  Millstone 

Hill  radar  measurements.  The  model  has  been  run  with  both  climatological 
neutral  winds  and  with  measured  winds.  The  model  results  with  adjusted 
winds  also  show  significant  improvement  in  Hmax  at  j^ht.  The  enhanced 
equatorw^d  winds  were  able  to  maintain  the  nighttime  ionosphere  (Nmax) 
under  active  magnetic  conditions  at  this  subauroral  site. 
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neutral  wind.  In  an  operational  sense,  there  will  be  no  additional  information  available  on 
the  ambient  levels  of  drifts  and  winds,  and  thus,  any  proposed  method  must  be  made  robust 
with  respect  to  independent  variations  in  both  quantities.  It  became  apparent  that  in  many 
cases  a  single  relationship  between  the  desired  ionospheric  parameter  existed  for  all  levels 
of  wind  and  drift,  while  in  other  cases  separate  relationships  were  defined  for  each  wind  or 
drift  level.  And  in  the  latter  case  a  separate  means  of  estimating  the  wind  or  drift  will  be 
required.  The  procedure  for  establishing  a  set  of  algorithms  for  deriving  Nmax  and  Hmax 
(and  secondarily,  effective  wind  and  effective  drift)  was  interactive,  and  was  performed  as 
follows: 

-  for  a  given  season,  level  of  solar  actiAuty,  and  local  time,  plot  Nmax  (or  Hmax)  as  a 
function  of  Ne(840)  at  a  range  of  latitudes  (each  plot  covering  all  wind  and  drift  levels  for 
one  latitude) 

-  flag  those  latitudes  where  no  knowledge  of  wind  or  drift  level  is  required 

-  at  other  latitudes,  determine  whether  separate  relationships  for  either  wind  or  drift 
level  (a  choice  made  from  visual  inspection) 

-  determine  a  latitude  where  an  estimate  of  wind  or  drift  can  be  made  purely  from  the 
value  of  Ne(840)  (specifically,  at  higher  latitudes,  vertical  drifts  have  little  influence,  and 
thus  estimates  of  wind  can  be  made-  similarly  at  the  equator,  the  vertical  drifts  dominate 
and  thus,  Ne(840)  usually  provides  a  guide  to  the  effective  drift  level). 

Such  a  set  of  coefficients  has  been  developed.  An  error  analysis  has  shown  that  the  standard 
deviation  of  the  algorithms  is  around  5-20%  in  Nmax,  and  slightly  less  in  Hmax.  Ignoring 
outliers,  the  standard  deviations  are  5-10%  in  both  Nmax  and  Hmax.  The  outliers  arise 
from  poor  estimates  made  for  effective  wind  or  drift,  and  thus,  it  is  preferable  to  base 
estimates  of  these  quantities  on  Ne(840)  values  from  more  than  one  latitude  in  each  case 
(assuming  that  the  wind  or  drift  level  applies  to  the  whole  region  and  not  just  locally,  that 
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the  variations  axe  physically  driven). 


4  Modeling  Applications  of  the  GTIM 

The  GTIM  determines  a  converged  solution  for  the  diurnal  variation  of  O'*"  density 
under  conditions  that  are  essentially  repeating  each  day.  Attempts  have  been  made  here 
to  adapt  the  model  for  some  level  of  real-time  responsiveness,  i.e.,  where  conditions  are 
changing  from  one  time  step  to  the  next.  Two  studies  of  this  type  have  been  conducted. 

4.1  Modeling  Impulsive  Disturbances 

Initially,  the  model  is  run  to  achieve  converged  conditions  for  what  is  basically  the 
quiescent  ionosphere  prior  to  the  disturbance  onset.  At  the  desired  onset  time  of  the  storm, 
a  different  set  of  Ap  values  are  invoked  to  drive  the  neutral  atmospheric  wind  and  density 
models  (these  drive  the  effects  at  middle  latitudes)  and  a  different  (perturbed)  diurnal 
vertical  drift  pattern  is  used  (driving  the  effects  at  low  latitudes)  along  the  lines  described 
in  Section  3  of  this  report.  The  perturbed  winds  and  drifts  are  applied  for  a  preset  time, 
and  the  Ap  history  subsequently  allows  the  model  ionosphere  to  relax  to  the  pre-onset 
conditions. 

Some  interesting  results  were' obtained  at  low  latitudes.  For  example,  at  latitudes  -10, 
0  and  -j-lO,  a  positive  phase  is  predicted  by  the  model  in  the  evening  sector  for  onset  times 
early  in  the  day.  Slight  negative  phases  are  predicted  at  -20  and  -1-20  degrees.  However, 
when  storm  onset  occurs  later  in  the  day,  there  is  an  evening  negative  phase  predicted  at  -10 
and  -1-10,  and  a  slight  positive  phase  at  -20  and  +20.  This  latter  point  is  shown  in  Figure  9, 
where  the  five  rows  show  the  variation  of  Nmax  with  local  time  at  the  five  named  latitudes. 
Leftmost  in  each  plot  is  the  quiet-time  variations,  and  this  is  repeated  as  a  dashed  line  for 
the  subsequent  days.  The  solid  line  shows  the  variation  in  Nmax  following  storm  onset,  and 
positive  and  negative  storm  phases  being  clearly  apparent.  Overall,  the  modeled 
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Figure  9.  Time  variations  of  modeled  Nmax  values  before  and  after  the 
onset  of  a  ma^etic  disturbance,  at  five  selected  latitudes.  Leftmost  in  eax:h 
plot  is  the  quiet-time  variations,  and  this  is  repeated  as  a  dashed  line  for  the 
subsequent  days.  The  solid  line  shows  the  variation  in  Nmax  following  storm 
onset,  and  positive  and  negative  storm  phases  being  clearly  apparent. 


ionospheric  storm  response  depended  both  on  the  local  time  at  onset  and  the  latitude  imder 
consideration. 

To  drive  the  model  for  middle  latitude  case  studies,  initial  conditions  are  first  achieved 
by  setting  a  constant  Ap  value  of  3.  At  storm  onset,  the  actual  Ap  history  of  the  storm 
to  be  modeled  is  utilized,  and  also  updated  each  three  hours  of  the  modeling  run  for 
subsequent  days.  At  middle  latitudes,  it  was  shown  that  storm-related  thermospheric 
changes  can  drive  significant  ionospheric  response  in  the  model.  There  is  often  a  good 
level  of  agreement  with  observations  of  storm-time  ionospheric  parameters  provided  by  the 
storm-driven  GTIM.  Consistently,  the  most  difficult  “features”  to  describe  by  the  model 
are  the  rapid  excursions  that  the  climatological  thermosphere  will  not  follow,  and  those 
occurring  in  predawn  hours  (when  the  theoretical  model  is  not  always  able  to  maintain  the 
ionosphere,  and  the  baseline  description  is  not  therefore  as  reliable).  The  storm  response 
in  the  model  is  a  fairly  rapidly  decreasing  function  of  latitude.  Overall,  the  level  of  broad 
agreement  attained  does  strongly  suggest  that  thermospheric  influences  are  dominant  in 
the  middle  latitude  disturbed  ionosphere,  and  that  even  the  climatological  model  MSIS  can 
provide  a  plausible  foundation  for  this  work. 

4.2  Modeling  Wavelike  Perturbations 

A  version  of  the  Phillips  Laboratory  GTIM  has  been  developed  that  utilizes  vertical 
drift  fields  with  wavelike  perturbations.  As  above,  the  model  is  first  allowed  to  converge 
to  quiescent  conditions,  before  at  a  predetermined  time  the  vertical  drifts  have  a  sinusoidal 
variation  imposed  upon  them. 

In  its  default  configuration,  the  GTIM  will  terminate  when  an  L-shell  reaches  a  preset 
minimum  value  (around  1.026).  This  needs  to  be  circumvented  if  one  is  to  study  waves  of 
long  periods  on  low  altitude  fieldlines.  In  essence,  the  method  currently  employed  is  one 
that  keeps  track  of  how  far  below  this  minimum  value  the  L-shell  is  apparently  drifting. 
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while  keeping  it  at  the  minimum  value  until  it  drifts  sufficiently  upward  to  compensate 
(equivalently,  there  is  a  lower  boundary  set  to  the  vertical  dependence  of  vertical  drifts,  with 
a  running  measure  taken  for  each  fieldline  to  ensure  that  no  net  drift  over  a  wave  period  is 
obtained).  A  figure  of  5%  is  considered  to  be  the  reliability  of  the  method  (determined  from 
numerical  tests  using  a  zero-amplitude  wave)  and  only  wave-related  effects  of  greater  than 
this  are  regarded  as  significant. 

Subsequently,  a  series  of  model  runs  have  been  conducted  appropriate  to  conditions 
where  a  16-day  periodicity  was  observed  in  ionospheric  data  (Forbes  and  Leveroni,  1992). 
Conditions  were  January,  solar  maximum  (1979),  US  sector,  and  in  the  data  some  20-30% 
(peak-peak)  variations  were  seen  in  noontime  Nmax.  The  series  of  model  runs  were  based 
on  several  different  onset-times  for  the  wavelike  variation  in  the  vertical  drift  pattern,  as 
well  as  for  two  different  amplitudes,  10%  and  25%  (20%  and  50%  peak-peak).  It  was  noted 
that  for  16-day  waves,  the  model  results  were  not  sensitive  to  the  onset-time  of  the  wave. 
Also,  the  amplitude  of  the  ionospheric  effects  were  not  linearly  related  to  the  amplitude  of 
the  input  disturbance. 

Figure  10  shows  one  means  of  displaying  the  results.  Here,  the  four  panels  show 
variations  in  Nmax  as  a  function  of  day  in  the  wave  for  selected  local  times  (0600LT, 

1200LT,  1800LT  and  2400LT).  Each  panel  shows  the  results  for  a  number  of  different 
(magnetic)  latitudes  (-18=solid,  -10=dotted,  0=dashed,  10=dot-dashed,  18=dot-dot-dot- 
dashed).  This  wave  had  a  peak-peak  amplitude  of  50%,  and  an  onset  time  of  1300LT. 

The  following  points  are  noted,  for  this  displayed  case,  and  in  general:  (1)  the 
corresponding  ionospheric  responses  are  largely  wavelike,  (2)  effects  are  largest  at  midnight, 
(3)  equatorial  effects  are  opposite  in  phase  to  AE  (equally,  AH),  (4)  at  local  times  where 
there  is  an  anomaly  crest,  effects  inside  and  outside  the  crests  are  opposite  in  sign  and  (5) 
effects  in  Hmax  are  much  smaller  than  in  Nmax. 
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Figure  10.  Wavelike  ionospheric  effects  in  the  GTIM  when  a  wavelike  E 
X  B  disturbance  is  applied.  This  case  was  made  for  January,  solair  maximum, 
Huancayo  sector  conditions,  using  a  16-day  wave  and  50%  peak-to-peak 
amplitude.  The  relative  variation  in  Nmax  at  selected  latitudes  (18=solid 
line,-10=dotted,0=dashed,-l-10=dot-dashed,-|-18=dot-dot-dot-dashed)  is 
plotted  as  a  function  of  the  day  during  the  period. 


Overall,  the  model  appears  to  give  a  good  qualitative  description  of  16-day 
periodicities  in  foF2,  if  one  creates  a  wavelike  disturbance  in  vertical  drift  fields. 
Quantitatively,  the  modeled  variations  are  typically  somewhat  smaller  than  those  observed. 
Better  quantitative  agreement  could  be  achieved  by  using  vertical  drifts  in  the  model 
that  were  obtained  from  actual  measurements,  in  order  to  better  describe  the  quiescent 
conditions. 

5  Applications  Studies 

Parallel  to  the  development  and  refinement  of  theoretical  ionospheric  models  is 
the  ongoing  development  of  applications-oriented  models  and  empirical  morphological 
descriptions  of  ionospheric  phenomena  and  their  terrestrial  impact.  Studies  of  this  latter 
type  are  addressing  particular  needs  of  users  of  real-time  or  monitoring  systems.  Several 
studies  of  this  type  are  described  in  the  following  section. 

5.1  Ionospheric  Ray- Tracing 

A  simple  profile  model  was  recently  developed  (Fox,  1994)  where  the  profiles  derived 
in  a  theoretical  model  at  low  latitudes  were  parametrized  in  terms  of  Chapman  functions 
for  each  of  the  topside  and  bottomside;  Chapman  representations  in  this  context  are  not  so 
much  physical  descriptions  as  convenient  mathematical  descriptions  of  a  basic  ionospheric 
layer  shape.  This  profile  model  differs  in  two  important  ways  (1)  the  scale  heights  of  the 
topside  and  bottomside  layers  are  allowed  to  be  height-dependent  (2)  the  profile  fit  is 
obtained  from  a  fit  to  both  peak  parameters  and  an  additional  control  parameter,  such  as 
slab  thickness  or  TEC.  The  model  can  thus  provide  climatological  profile  estimates,  as  well 
being  able  to  accommodate  additional  real-time  measurements. 

The  task  undertaken  here  was  to  adapt  this  profile  model  for  raytracing  applications. 
The  general  problem  in  radiowave  propagation  is  to  be  able  to  predict  how  various 
frequencies  of  waves  will  propagate  for  a  given  circuit  under  given  conditions.  Systems 
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already  exist  that  can  provide  estimates  based  on  climatology  but  it  will  not  be  clear  in  the 
field  how  any  additional  (real-time)  observations  may  be  used  to  help.  The  model  described 
above  is  thus  potentially  very  useful. 

In  the  first  instance,  we  used  the  model  to  provide  climatological  profiles  for  a 
number  of  locations  along  a  designated  circuit,  and  used  a  two-dimensional  ray-tracing 
model  to  define  the  (zeroeth  order)  propagation.  Note  that  in  this  profile  model  both 
the  electron  density  profiles  and  their  first  derivatives  are  analytic.  Next,  over  the  array 
of  frequencies  considered,  the  elevation  angle  required  to  reach  the  receiver  from  the 
transmitter  is  determined.  Plots  of  elevation  angle  versus  frequency  are  used  as  the  basis 
for  comparisons.  Next,  a  ground  truth  ionosphere  was  set  up  for  the  same  conditons,  and 
the  profile  model  given  different  amounts  of  information  (peak  or  topside  information,  at 
various  sets  of  locations)  to  update  the  specification  and  to  ray-trace.  What  is  being  tested 
here  is  the  ability  of  the  Semi-empirical  Improved  Middle  Plus  Low-latitude  Ionosphere 
(hereafter,  SIMPLI)  model  to  assimilate  a  small  amount  of  data  about  the  circuit  and  to 
successfully  update  radiowave  propagation  predictions  accordingly.  The  means  of  testing 
is  to  compare  the  elevation  angle  loci  of  different  ionospheres  where  each  has  resulted 
from  a  SIMPLI  description  using  different  amounts  of  “real”  (in  this  case,  from  the  Bent 
[1973]  model)  data.  Imagine  first  only  the  Nmax  value  at  the  circuit  midpoint  is  known; 
different  methods  of  updating  the  description  of  profile  parameters  along  the  raypath  were 
tested.  Successively,  the  same  tests  were  performed  when  other  parameters  (Hmax  and 
slab  thickness)  and/or  specified  at  other  locations  (along  the  circuit)  are  known.  This 
determined  how  much  data  is  required  for  the  SIMPLI  ray-tracing  model  to  adequately 
revise  its  propagation  predictions. 

We  illustrate  typical  results  with  the  circuit  Tx=(30,30),  Rx=(20,20),  June, 

FIO. 7=150,  1200UT,  presented  in  Figure  11.  The  solid  line  is  the  initial  guess  from  SIMPLI 
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SIMPLI  raytracing  :  Nmax  only 


Figure  11.  A  demonstration  of  the  abilities  of  the  SIMPLI  profile  model 
in  radiowave  propagation  predictions.  Here,  the  solid  line  indicates  the 
elevation  angle  versus  transmission  frequency  for  a  given  circuit  obtained 
from  climatological  estimates  (electron  density  contours  in  the  left  insert). 
The  dotted  line  shows  the  same  predictions  made  from  the  “ground-truth” 
ionosphere,  shown  in  the  right  insert.  The  three  topmost  lines  show  the 
results  when  Nmax  is  updated  at  the  circuit  mid-point,  end-points,  and  both, 
and  when  SIMPLI  parameters  along  the  path  axe  shifted  only  the  zeroeth 
order.  The  three  lines  lying  nearest  to  the  dotted  line  correspond  to  Nmax 
being  refit  along  the  whole  path. 
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for  the  locus  of  elevation  angle  versus  frequency  for  the  circuit,  based  on  the  ionosphere 
seen  in  the  left  contour  plot  included  in  the  figure.  The  dashed  line  is  the  locus  from  the 
“actual”  ionosphere  (from  the  Bent  model),  seen  in  the  right  contour  plot,  and  clearly  very 
different  to  the  SIMPLI  climatology  in  this  case.  Note  both  the  different  elevation  angles  in 
the  true  case,  as  well  as  the  greater  range  of  usable  frequencies.  It  is  important  to  be  able 
to  track  changes  of  this  nature.  Also  seen  in  this  figure  are  two  sets  of  three  different  lines, 
where  the  thinner  lines  are  derived  from  one  method  of  updating  profile  parameters,  and  the 
thicker  lines  from  another,  and  where  the  data  used  is  Nmax  from  the  midpoint  (dashed), 
ends  of  the  circuit  (dot-dashed)  and  all  of  the  above  (dot-dot-dot-dashed).  In  this  case,  the 
“fitting”  method  is  seen  to  give  a  fairly  good  update  for  the  ray  propagation,  even  when 
only  Nmax  at  the  midpoint  is  known. 

This  same  type  of  test  was  peformed  over  a  variety  of  one-hop  circuits  over  a  variety  of 
conditions.  In  general,  the  fitting  technique  provided  the  more  reliable  means  of  updating 
the  SIMPLI  model  for  ray-tracing  applications  and  overall,  good  results  can  be  achieved 
when  only  Nmax  at  the  circuit  midpoint  is  available  as  an  update.  Nmax  available  at  either 
the  transmitter  or  receiver  location  also  provides  sufficient  information  to  substantially 
improve  propagation  predictions. 

I 

5.2  Sensing  the  F  region  fifoin  DMSP 

Earlier  work  performed  at  Boston  University  established  a  means  of  inferring  Nmax 
and  Hmax  values  from  (GTIM  modeled)  observations  of  nightglow  intensities  (of  the  type 
to  be  measured  by  the  Special  Sensor  Ultraviolet  Spectroscopic  Imager,  or  SSUSI).  Here, 
a  further  study  has  been  made,  involving  the  use  of  these  Nmax  and  Hmax  estimates  in 
conjunction  with  direct  measurements  of  the  in  situ  electron  density  from  the  Special  Sensor 
for  Ion,  Electron  and  Scintillation  Densities  (SSIES).  The  specific  task  tackled  was  to  derive 
a  means  of  combining  the  in  situ  and  peak  ionospheric  parameters  into  a  specification  of 
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the  entire  profile  (or  the  total  electron  content,  or  TEC),  performing  data  fusion  on  the 
ionospheric  inputs  from  different  sensors. 

The  method  adopted  is  briefly  as  follows  and  is  based  on  the  profile  model  described  in 
the  previous  section  of  this  report.  Chapman-like  layers  of  altitude-dependent  scale  heights 
are  anchored  to  the  peak  parameters.  Initial  estimates  of  the  scale  heights  axe  obtained 
from  a  parametrization  of  results  from  the  GTIM.  The  resultant  in  situ  (840  km)  density  is 
compared  to  the  SSIES  value,  and  the  difference  is  used  to  alter  the  internal  scale  heights, 
until  a  converged  profile  solution  is  obtained.  In  this  fashion,  smoothly-varying  profiles 
are  derived  that  axe  anchored  exactly  at  the  peak  and  to  a  predetermined  precision  in  the 
topside.  Estimates  of  TEC  made  thusly  should  be  reliable  because  the  profile  will  always 
asymptote  to  zero  at  high  altitudes.  It  should  also  be  noted  that  although  the  profile  model 
is  based  on  there  being  one  topside  “layer”,  it  is  not  limited  to  applications  where  only  one 
ionic  species  is  present.  The  ability  to  adjust  the  effective  scale  height  compensates  for  the 
transition  from  O"’’  to 

Documentation  has  been  prepared  for  this  procedure  and  has  been  delivered  to  the 
Applied  Physics  Laboratory  of  Johns  Hopkins  University,  in  Laurel,  MD. 

5.3  Scintillations 

I 

A  review  was  prepared  on  the  scintillation  activity  reported  in  the  literature  at  GPS 
frequencies,  including  observations  at  1.2  GHz  to  7  GHz  with  the  appropriate  scaling  for 
the  higher  frequencies.  It  also  included  available  data  from  GPS  scintillation  measurements. 
A  first  summary  of  available  data  was  given  at  the  Institute  of  Navigation-GPS  Meeting  in 
Salt  Lake  City  in  the  Spring  of  1995;  this  report  has  been  published  in  the  Proceedings  of 
the  Institute. 

The  most  severe  scintillation  occurs  in  the  equatorial  anomaly  region.  In  this  region, 
during  the  solar  maximum  periods,  amplitude  scintillations  at  1.5  GHz  may  exceed  20  dB 
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for  several  hours  after  sunset.  Ascension  Island  in  the  Atlantic,  Diego  Garcia  in  the  Indian 
Ocean,  Hong  Kong  and  Taiwan  in  the  Pacific  are  some  of  the  stations  that  fall  directly 
under  the  anomaly  region.  However,  areas  within  a  few  degrees  of  the  magnetic  equator 
may  show  only  5-7  dB  fades.  The  occurrence  of  scintillation  in  the  equatorial  region  is  also  a 
function  of  longitude. 

The  other  potentially  active  regions  are  at  auroral  and  polar  cap  latitudes.  In  the 
central  polar  cap  in  years  of  solar  maximum,  GPS  receivers  may  suffer  >10  dB  fades. 

During  magnetic  storms  auroral  latitude  reception  may  experience  scintillation  at  a  level 
to  produce  loss  of  lock.  Models  of  scintillation  occurrence  exist  and  they  are  of  use  for 
planning  purposes.  However,  the  day-to-day  variations  of  ionospheric  parameters  malce 
short  term  forecasting  difficult. 

GPS  satellites  offer  a  unique  source  for  measurements  of  amplitude  and  phase 
scintillation  on  a  global  scale.  One  receiver  can  record  scintillation  magnitudes  and 
spectra  at  multiple  propagation  paths  in  the  overhead  sky.  The  data  can  be  used  to  study 
ionospheric  plasma  structures,  develop  weather  models  of  scintillation  and  can  be  scaled 
in  frequency  to  support  many  operational  systems.  In  turn,  real  time  observations  could 
allow  for  warning  of  fading  problems  since  a  laxge  portion  of  the  sky  would  be  covered 
in  the  reception  of  signals  from  the  multitude  of  satellites  in  the  system.  An  approach 
to  a  warning  system  has  been  developed.  Any  interpretation  of  the  fading  structure  and 
scintillation  spectra  must  take  into  account  both  the  ionospheric  motion  and  the  velocity 
of  the  satellite.  For  GPS  receivers  in  aircraft,  the  motion  of  the  heading  of  the  aircraft  also 
needs  to  be  considered.  With  the  use  of  more  data  and  an  understanding  of  the  behavior  of 
equatorial  irregularities,  techniques  could  be  used  to  minimize  the  problem. 
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6  Meridional  Wave  Structures 


SETA  density  measurements  (normalized  to  a  nominal  altitude  of  200  km)  were 
analyzed  to  delineate  characteristics  of  relatively  small-scale  waves  (horizontal  structures  of 
order  150  to  2500  km)  in  the  data.  (Note:  This  analysis  is  restricted  to  in-track  (latitudinal) 
structures;  the  analysis  of  longitudinal  structures  is  addressed  further  below.)  The  spectral 
analyses  were  applied  to  all  available  orbits  (approximately  15  to  18  per  day)  during  6  days 
in  1983  and  10  days  in  1982,  covering  a  range  of  geomagnetic  conditions. 

The  complete  results  of  this  project  appear  in  Forbes  et  al.  (1995).  A  brief  summary 
of  results  is  presented  below,  and  descriptions  of  the  equations  and  algorithms  (with 
samples  of  results)  used  in  the  analysis  are  provided  in  the  following  sections. 

6,1  Summary  of  Results 

As  a  result  of  this  analysis,  the  lower  thermosphere  is  found  to  be  frequently 
characterized  by  wavelike  structures  corresponding  to  peak-to-peak  density  variations  of 
order  10-40%.  Orbit  by  orbit  spectral  analyses  of  densities  normalized  to  200  km  altitude 
were  performed,  and  these  are  combined  together  to  construct  average  periodograms  for 
magnetically  quiet  and  active  periods.  Integration  of  these  curves  provides  the  relative 
percentage  of  total  density  variance  attributable  to  various  horizontal  scales;  for  instance, 
10%  for  Aif  <  750  km  and  70%  for  1500  <  \h  <  2500  km.  As  geomagnetic  activity 
increases  from  quiet  (Kp  <  2)  to  active  {Kp  >  4)  levels,  the  spectral  energy  increases 
by  a  factor  of  two  for  wavelengths  in  the  range  150  -  500  km  and  by  a  factor  of  four  at 
longer  wavelengths  (1000  -  2500  km).  In  part  the  spectral  energy  enhancement  at  the 
longer  wavelengths  reflects  mesoscale  density  variations  unrelated  to  wave  activity,  but 
“direct”  waves  generated  by  auroral  activity  and  which  are  known  to  travel  long  distances 
probably  also  make  important  contributions  at  these  wavelengths.  Similarly,  a  comparison 
of  day  vs.  night  periodograms  for  all  orbits  during  1983  was  performed,  illustrating  the 
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enhancement  in  wave  activity  which  occurs  at  scales  shorter  than  400  km  during  the  night, 
whereas  the  spectral  densities  at  wavelengths  longer  than  500  km  are  virtually  the  same. 
The  magnetic  activity  enhancement  at  shorter  scales  is  also  significant,  especially  at  low 
latitudes,  since  waves  with  these  horizontal  wavelengths  are  not  expected  on  theoretical 
grounds  to  propagate  very  far  within  the  thermosphere  away  from  an  auroral  source. 
Moreover,  the  shorter  scale  waves  regularly  occur  during  geomagnetically  quiet  times  when 
an  auroral  source  mechanism  is  unlikely. 

Much  of  our  efforts  were  devoted  to  developing  and  applying  various  spectral  methods 
to  the  data  set.  The  following  describes  the  methods  of  analysis  which  were  developed,  and 
presents  sample  results  from  the  analysis  which  are  alluded  to  in  the  above  summary. 

6.2  Method  of  Spectral  Analysis 

Standard  periodogram  spectral  analyses  were  applied  to  each  orbit  of  uninterrupted 
data.  The  orbits  were  categorized  into  three  levels  of  magnetic  activity  in  order  to 
determine  the  dependence  of  periodicities  on  the  magnetic  activity.  The  categories  are  (a) 
quiet;  Kp  <  2;  (b)  average:  2  <  Kp  <  4;  (c)  active:  Kp  >  4.  In  this  section,  we  describe  the 
results  of  periodogram  spectral  analysis  performed  on  the  density  data. 

The  periodogram  spectral  analysis  is  based  on  estimating  the  distribution  in  frequency 

1 

of  the  power  Px{iv)  or  the  power  spectral  density  (PSD)  of  a  random  signal  xfn]  (density 
measurements  in  this  analysis).  This  quantity  can  be  shown  to  be  related  to  the  squared 
magnitude  of  the  Fourier  transform  of  the  signal  x{n]  by  an  extension  of  the  Wiener- 
Khinchin  theorem  as  follows: 

P.(to)  =  E  |e" -«  xln]e->«"f 

where  E[.]  denotes  the  expectation  operator.  Neglecting  the  expectation  operator,  the  PSD 
of  the  signal  x[n]  of  length  N  may  be  estimated  as: 

hiw)  =  jr  eS 
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which  is  the  squared  magnitude  of  the  Fourier  transform  of  x\n]  and  can  be  evaluated 
using  the  FFT  algorithm. 

PSD’s  were  generated  for  each  uninterrupted  orbit  in  the  data  base  using  4096-point 
FFT’s.  Since  each  orbit  of  the  satelite  takes  approximately  90  minutes,  1296  data  points 
were  analysed  for  each  orbit  in  the  case  of  the  1983  measurements. 

6.3  Periodogram  Averaging 

The  power  spectra  can  be  estimated  more  “reliably”  by  attempting  to  approximate 
the  expectation  operator  described  in  the  previous  section.  The  technique  is  known  as 
periodogram  averaging  and  it  is  described  in  the  following. 

Suppose  that  ^i[n],  X2Wj  K  uncorrelated  instances  (data  records)  of 

the  same  random  signal  x[n].  The  averaged  periodogram  is  then  defined  as:  Pavg{w)  = 

^  Z)m=i  Pm{w) 

A  I  *12 

where  Pm  is  the  m**  periodogram  given  by:  Pmiw)  =  ^ 

which  is  the  squared  magnitude  of  the  Fourier  transform  of  the  data  record 
The  variance  of  the  averaged  periodogram  can  be  shown  to  be  decreased  by  a  factor  of 
K  while  its  expected  value  remains  the  same  as  the  true  PSD.  In  this  sense  the  averaged 
periodogram  is  a  more  reliable  estimate  of  the  power  spectrum. 

Our  averaged  periodogram  analyses  are  twofold.  First,  the  average  of  the  spectra  for 
all  of  the  uninterrupted  orbits  in  each  day  is  calculated.  Second,  the  orbits  are  categorized 
into  three  levels  of  magnetic  activity  as  described  in  the  previous  section,  and  for  each  level 
of  magnetic  activity  an  averaged  periodogram  is  calculated. 

The  averaged  periodogram  corresponding  to  quiet  and  active  orbits  during  the  six 
days  in  1983  are  illustrated  in  Figure  12.  This  illustrates  the  enhancement  in  wave  activity, 
especially  at  wavelengths  longer  than  1,000  km,  which  is  generated  in  conjunction  with 
magnetic  disturbances. 
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Similarly,  a  comparison  of  day  vs.  night  periodograms  for  all  orbits  during  1983 
are  illustrated  in  Figure  13.  This  plot  illustrates  the  enhancement  in  wave  activity  which 
occurs  at  scales  shorter  than  400  km  during  the  night,  whereas  the  spectral  densities  at 
wavelengths  longer  than  500  km  are  virtually  the  same. 

6.4  Wavelength  Contribution  to  Signal  Energy 

In  order  to  gain  quantitative  insight  into  the  contribution  of  various  wavelength  ranges 
to  the  total  signal  energy,  an  equation  known  as  Parseval’s  relation  may  be  utilized,  which 
relates  a  signal’s  total  energy  to  its  total  spectral  energy  by: 

E^-oo  kWP  =  \X{w)\^dw 

where  a:[n]  is  the  discrete-time  signal  corresponding  to  the  data  values,  and  X{w)  is  the 
spectrum  of  the  signal.  The  quantity  on  the  left-hand  side  of  the  equation  is  the  total 
energy  in  the  time-domain  signal  x[n],  and  |X(ta)p  is  referred  to  as  the  energy-density 
spectrum.  The  quantity  of  interest  here,  which  we  refer  to  as  the  “percentage  contribution” , 
is  the  ratio  of  the  energy  in  an  arbitrary  wavelength  range  to  the  total  signal  energy. 
According  to  Parseval’s  relation  this  quantity  is  proportional  to  the  ratio  of  the  energy 
in  the  arbitrary  wavelength  range  to  the  total  spectral  energy.  Since  the  wavelengths  of 
interest  here  are  all  in  the  range  of  0  to  2500  km,  the  “total”  energy  may  be  replaced  by  the 
energy  in  the  range  of  0  to  2500  km.  The  percentage  contribution  may  therefore  be  defined 
as  the  ratio  of  the  energy  in  an  arbitrary  wavelength  range  to  the  total  energy  of  the  signal 
with  wavelengths  between  0  and  2500  km.  The  averaged  periodogram  results  described  in 
the  previous  section  are  used  to  compute  the  percentage  contribution  for  the  density  data 
during  1983  as  a  continuous  measure  of  the  ratio  for  all  wavelengths  in  the  range  of  0  to 
2500  km.  The  percentage  contributions  of  the  wavelengths  are  therefore  both  equal  to  0% 
at  0  km  and  100%  at  2500  km,  and  increase  continuously  with  increasing  wavelength.  As 
noted  previously,  integration  of  these  curves  provides  the  relative  percentage  of  total  density 
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variance  attributable  to  various  horizontal  scales;  for  instance,  10%  for  \ji  <  750  km  and 
70%  for  1500  <  Aff  <  2500  km.  As  geomagnetic  activity  increases  from  quiet  (Kp  <  2)  to 
active  {Kp  >  4)  levels,  the  spectral  energy  increases  by  a  factor  of  two  for  wavelengths  in  the 
range  150  -  500  km  and  by  a  factor  of  four  at  longer  wavelengths  (1000  -  2500  km). 

The  following  conclusions  may  be  drawn  from  the  above  work:  Meridional  wave 
structures  exist  in  thermospheric  density  with  specific  characteristics  which 
depend  on  local  time  and  magnetic  activity.  These  variations  are  not  included  in 
empirical  models.  The  majority  of  density  variance  associated  with  these  waves 
occurs  at  scales  in  excess  of  1,000  km,  and  therefore  in  principle  their  presence 
can  be  replicated  in  first-principles  models  of  the  thermosphere-ionosphere 
system. 

7  Magnetic  Storm  Perturbations 

SETA-2  and  SETA-3  data  have  been  analyzed  during  the  following  periods  containing 
significant  changes  in  geomagnetic  activity: 

(a)  July  20-27,  1983 

(b)  August  2-10,  1983 

(c)  August  2-9, 1982 

(d)  September  3-11,  1982 

(e)  December  4-14,  1983 

The  results  of  this  study  are  reported  fully  in  Forbes  et  al.  (1996).  A  summary  of  results  is 
presented  below,  followed  by  sample  results  for  two  storm  periods  which  are  characteristic  of 
the  others  studied. 
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7.1  Summary  of  Results 

Density  measurements  near  200  km  from  the  SETA  (Satellite  Electrostatic  Triaxial 
Accelerometer)  experiment  are  used  to  delineate  the  temporal,  seasonal-latitudinal  and  day 
vs.  night  dependences  of  the  response  to  magnetic  storm-related  perturbations  in  high- 
latitude  energetic  inputs.  Five  periods  of  geomagnetic  activity  are  analyzed,  and  yield 
consistent  results  which  can  be  interpreted  within  the  framework  of  recent  thermosphere- 
ionosphere  simulations.  In  response  to  a  change  in  magnetic  activity  level  from  quiet  {Kp 
«  1-2)  to  active  [Kp  k  4-7)  conditions,  an  increase  in  daytime  (1030  LT)  density  of  order 
50-70%  occurs  between  60  and  80°  geographic  latitude  in  the  summer  hemisphere,  with 
about  half  the  maximum  response  in  the  winter  hemisphere.  This  difference  is  mainly  due 
to  different  ionization/conductivity  levels  (and  hence  joule  heating  rates)  between  the 
hemispheres.  Dayside  penetration  of  the  disturbance  at  about  the  50%  intensity  level  is 
realized  at  the  equator,  whereas  in  the  winter  hemisphere  equatorward  penetration  is  much 
weaker.  These  effects  are  connected  with  the  prevailing  solar-driven  circulation;  the  net 
summer-to- winter  meridional  flow  facilitates  equatorward  advection  of  the  disturbance 
bulge  in  the  summer  hemisphere,  but  hinders  advection  in  the  winter  hemisphere.  In  both 
hemispheres  the  daytime  solax-driven  diurnal  circulation  tends  to  oppose  equatorward 
penetration  to  the  same  degree.  However,  on  the  nightside  (2230  LT)  penetration  at  nearly 
the  100%  level  of  both  summer  and  winter  disturbance  bulges  are  realized  to  within  20° 
of  the  geographic  equator.  This  behavior  is  associated  with  the  equatorward  advection  in 
both  hemispheres  consistent  with  the  nighttime  component  of  the  solar-driven  circulation. 
Comparisons  with  the  MSISE90  model  \Hedin,  1991]  show  it  to  capture  the  salient  features 
of  the  daytime  behavior,  but  exhibits  little  day /night  asymmetry,  in  contrast  to  the 
experimental  results. 
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7.2  Sample  Result 


Various  data  visualization  tools  have  been  developed  to  depict  the  temporal,  day 
vs.  night,  and  latitudinal/seasonal  trends  characterizing  the  response  of  thermospheric 
mass  density  to  changes  in  high-latitude  heating  during  magnetically  active  periods.  Data 
between  170  and  240  km  are  normalized  to  a  constant  altitude  of  200  km.  Sliding  11-rev 
(about  16  hours)  means  are  applied  to  the  data  to  remove  smaller-scale  features,  and  to 
serve  as  a  “trend”  for  delineation  of  the  residual  densities  (see  following  section). 

Kp  variations  for  two  of  the  study  periods,  July  20-27,  1983,  and  August  1-10,  1983, 
are  depicted  in  Figure  14.  Figure  15  illustrates  the  latitude-time  variations  in  total  mass 
density  “trend”  occurring  during  the  July  1983  period.  Note  the  much  larger  response 
occurring  in  the  summer  hemisphere;  this  is  due  in  part  to  the  enhanced  Joule  heating 
connected  with  the  larger  ionization  levels  characterizing  the  summer  ionosphere.  There 
is  also  the  tendency  for  the  density  perturbation  to  be  carried  from  summer  to  winter  due 
to  the  prevailing  thermospheric  circulation  in  that  direction.  The  combination  of  these  two 
effects  gives  rise  to  the  large  seasonal  asymmetry  in  the  density  response  depicted  in  Figure 
15. 

Figure  16  illustrates  the  SETA  density  response  in  the  nightside  (2230  LST)  for 
the  July  1983  period.  Note  that  the  density  response  penetrates  much  more  easily  to  the 
equatorial  region;  this  is  due  to  the  enhanced  equatorward  winds  which  characterize  the 
nighttime  thermosphere  as  opposed  to  the  daytime  thermosphere.  The  corresponding  figure 
for  MSISE90  (not  shown)  does  not  reflect  this  behavior,  pointing  to  a  shortcoming  of  this 
empirical  model. 

The  following  conclusions  may  be  drawn:  The  response  of  thermospheric  density 
to  magnetic  activity  is  characterized  by  dependences  on  local  time  and  season 
which  may  be  interpreted  in  terms  of  prevailing  (solar-driven)  circulation 
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Figure  15.  Daytime  density  variations  during  July  20-27,  1983. 


systems.  The  data  exhibit  several  characteristics  which  are  not  taken  into 
account  in  existing  empirical  models.  It  is  recommended  that  additional  storms 
be  analyzed  to  further  solidify  these  results,  and  that  empirical  models  used  by 
the  Air  Force  be  updated  accordingly. 

8  Longitudinal  Wave  Structures 

8.1  Longitude/UT  Variation 

The  existence  of  a  “longitude/UT”  effect  in  the  MSIS  series  of  empirical  models, 
which  numerically  describes  the  geomagnetic  control  of  high-latitude  thermospheric  density, 
suggests  that  we  examine  the  SETA  data  for  such  variations.  Since  the  MSIS  models  suffer 
from  sparse  data  coverage  at  low  orbital  altitudes  and  high  solar  activity,  this  comparison 
should  be  valuable  in  ascertaining  the  inadequacies  of  empirical  models  under  these 
conditions.  This  comparison  has  been  a  major  focus  of  this  year’s  contractual  activities. 

Figure  18  illustrates  the  time  variation  of  SETA  densities,  normalized  to  200  km  and 
averaged  between  20  and  40  degrees  latitude,  during  July  21-26,  1983.  MSIS90  density 
values  were  calculated  for  each  SETA  data  point,  and  averaged  in  an  identical  fashion; 
these  values  are  shown  by  the  dotted  line.  Note  that  both  the  data  and  model  exhibit 
a  repeatable  cycle  with  period  of  24  hours  or  360  degrees  longitude;  it  is  this  variation 
which  we  refer  to  as  the  longitude/UT  effect.  Note  that  there  is  the  appearance  of  the 
longitude/UT  “amplitude”  increasing  with  magnetic  activity  (a  magnetic  storm  occured  on 
July  24). 

For  the  present  study,  we  primarily  confine  ourselves  to  the  magnetically  quiet  periods 
during  July  20-24,  1983,  and  December  4-10,  1983.  All  comparisons  between  MSIS  and 
SETA  correspond  to  SETA  densities  normalized  to  200  km  with  the  longitude/UT  effect 
removed  from  the  MSIS  model.  Variations  in  SETA  data  which  we  compare  with  MSISE90 
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Fi^re  18.  Compajison  of  SETA  and  MSISE90  density  variations 
averaged  over  20-20  degrees  latitude  at  1030  LT  during  July  20-26, 1983. 


were  extracted  as  follows.  First,  the  long-term  trend  was  established  by  performing  a  17- 
orbit  running  mean.  The  residuals  from  this  mean  were  then  smoothed  (5-orbit  running 
mean)  to  remove  the  short-term  variations.  The  result  of  this  process  is  indicated  in  Figure 
19  for  dayside  (1030  LT)  data  averaged  over  the  0-10°  latitude  range.  The  solid  line  in 
Figure  19  represents  the  least-squares  fit  of  a  sinusoid  with  a  period  of  24  hours  (360° 
in  longitude)  to  the  SETA  density  residuals  (circles).  The  dashed  curve  represents  the 
corresponding  MSISE90  prediction  for  this  period.  Note  that  the  amplitude  is  given  in 
terms  of  percent  variation  about  the  mean  density.  Figure  19  reveals  a  distinct  periodicity 
in  both  the  SETA  data  and  MSIS  which  we  are  referring  to  here  as  the  “longitude/UT” 
effect,  although  longitude  and  UT  variations  cannot  be  separated  or  distinguished  from  the 
sun-synchronous  SETA  orbit.  The  MSIS  model  exhibits  a  slightly  smaller  ampltiude  and  a 
slight  shift  in  phase  with  respect  to  the  SETA  density  oscillation. 

Figure  20  summarizes  results  from  the  above  procedure  for  all  available  latitudes  for 
the  daytime  data  during  this  period.  Amplitudes  are  calculated  (a)  taking  into  account  the 
magnetically-active  period,  and  (b)  considering  only  the  quiet  days.  The  major  model/data 
discrepancy  to  be  gleaned  from  this  plot  is  the  signficiant  latitudinal  asymmetry  depicted  by 
the  SETA  data  -  a  factor  of  2  to  3  enhancement  of  the  oscillation  in  the  S.  Hemisphere;  on 
the  other  hand  the  MSIS  model  is  symmetric  about  the  equator. 

An  alternative  form  of  depiction  is  presented  in  Figures  21  and  22.  Here,  we  provide 
grey-scale  plots  of  the  density  comparisons  which  also  take  into  account  the  phase 
(i.e.,  longitude  of  maximum)  of  the  “longitude/UT”  oscillation.  Only  the  magnetically 
undisturbed  data  are  used  in  the  fits.  These  plots  also  provide  comparisons  for  the 
nighttime  data  for  the  July  period,  and  daytime  data  for  the  December  period.  For  the 
daytime  data  during  July  (left-hand  panels  in  Figure  21),  SETA  and  MSIS  compare  well  in 
phase  (i.e.,  longitude  of  maximum),  but  the  SETA  data  suggest  the  phenomenon  to 
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Figure  19.  Comparison  of  ’’longitude/UT”  variation  in  SETA  densities 
(circles)  with  MSISE90  (dotted  line)  during  the  magnetically  quiet  period  of 
July  21-23,  1983.  Solid  line  represents  least-squares  fit  of  sinusoid  to  SETA 
data. 
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Figure  20.  Comparison  of  ampltiudes  of  ”longitude/UT”  variation 
between  SETA  and  MSISE90  as  a  function  of  latitude,  for  daytime  data 
during  July  21-23,  1983. 
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Figure  21.  Longitude/UT  dependence  revealed  by  SETA  (top)  and 
MSISE90  (bottom)  total  mass  densities  at  200  km  during  July  21-23,  1983. 
Left:  daytime  («1030  LT).  Right:  nighttime  («2230  LT). 
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Figure  22.  Longitude/UT  dependence  eis  revealed  in  SETA  (left)  and 
MSISE90  (right)  daytime  («1030  LT)  total  mass  densities  at  200  km  during 
December  5-10,  1983. 


be  much  more  robust  in  the  Southern  Hemisphere  as  compared  with  MSIS.  The  right-hand 
panels  illustrate  the  same  results,  except  for  nighttime  (approx  2230  LT).  Here,  the  converse 
is  true:  SETA  is  more  symmetric  about  the  equator,  and  MSIS  is  more  asymmetric,  with 
much  larger  perturbations  in  the  Northern  Hemisphere.  Figure  21  illustrates  similar 
comparisons  for  daytime  data  during  December  1983,  illustrating  a  very  interesting  180 
degree  shift  (for  both  MSIS  and  SETA)  in  the  phase  (longitude  of  maximum)  of  the 
phenomenon,  as  compared  to  the  July  period.  Physical  mechanisms  underlying  the  above 
behaviors,  including  the  seasonal  phase  shift,  are  not  understood  at  the  present  time. 

The  main  conclusion  to  be  drawn  from  these  comparisons  is:  There  exist  significant 
seasonal/latitudinal  differences,  in  amplitude  and  phase,  between  large-scale 
longitudinal  structures  revealed  by  SETA  data  in  comparison  to  the  MSISE90 
model  during  magnetically  undisturbed  periods;  these  differences  should 
be  incorporated  in  future  empirical  modelling  efforts  in  order  to  improve 
specifications  of  thermospheric  density  for  various  Air  Force  applications. 

8.2  Longitudinal  Wave  Structures 

Given  the  above  large-scale  trends  in  the  longitude  structure  of  thermospheric  density, 
subsequent  efforts  have  focused  on  analysis  and  interpretation  of  the  “residual  density 
variations”  superimposed  on  the  large-scale  structure.  Samples  of  these  residual  structures 
for  daytime  data  during  July  20-24,  1983  and  December  5-10,  1983,  are  provided  in  Figure 
23.  It  is  evident  that  there  exist  band-like  structures  extending  in  the  north-south  direction. 
The  following  describes  our  efforts  to  characterize  these  features  quantitatively. 

The  main  difficulty  faced  in  analyzing  and  interpreting  these  features  is  with  the 
space-time  sampling  inherent  from  a  satellite.  The  SETA  satellites  are  in  sun-synchronous 
orbits,  which  means  that  local  time  is  fixed  and  the  satellite  track  is  nearly  at  an  angle 
of  90°  with  respect  to  the  equator  (actually,  the  angle  is  about  83°).  The  earth  (and 
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atmosphere)  is  thus  rotating  below  the  satellite  at  an  approximate  rate  of  360°  longitude 
per  day.  Thus,  time  and  longitude  are  varying  continuously,  and  it  is  therefore  difficult 
to  separate  longitude  variations  from  time  variations  in  the  atmospheric  density.  This  is 
especially  the  case  when  wave  features  may  in  principle  be  propagating  eastward  and/or 
westward  with  resepct  to  the  surface  of  the  earth,  as  we  well  know  that  case  to  be  in 
the  troposphere  and  stratosphere.  Spectral  analysis  of  the  data  therefore  leads  to  some 
ambiguity  in  interpretation.  The  difficulty  may  be  summarized  through  the  relationship 
given  below 

fcs  =  m  ±  1/r 

where  ks  is  the  zonal  wavenumber  viewed  from  the  frame  of  the  satellite,  m  and  T  are  the 
earth-fixed  zonal  wavenumber  (integer)  and  period  (in  days)  for  a  particular  oscillation, 
and  the  negative  (positive)  sign  corresponds  to  westward  (eastward)  propagation,  kg  need 
not  be  an  integer.  So,  if  kg  is  found  to  be  1.33,  for  instance,  from  straightforward  spectral 
analysis  of  the  data,  this  could  be  interpreted  as  a  westward-propagating  m=2  oscillation 
with  T=1.5  days,  or  an  eastward  traveling  wave  with  m=l  and  T=3  days.  Other  higher- 
order  combinations  are  in  theory  possible,  and  often  it  is  necessary  to  either  live  with  the 
ambiguity,  or  to  deduce  the  physically  correct  answer  by  imposing  other  constraints. 

The  most  prominent  value  of  kg  derived  from  spectral  analyses  of  the  SETA  at 
different  latitudes  was  a  value  of  4,  with  secondary  peaks  often  occurring  at  2,  3, 1,  and  5. 

A  superposition  of  spectra  for  the  various  latitude  bands  in  the  N.  Hemisphere  for  the  July 
daytime  data  (cf.  Figure  23)  is  illustrated  in  Figure  24.  For  kg  =  2,  for  instance,  we  have  as 
the  lowest-order  mathematically  plausible  waves  (m,  T)  =  (0,  1/2),  (2,oo),  (1,-1).  The  first 
corresponds  to  a  zonally-symmetric  semidiurnal  tide;  the  second  to  a  quasi-stationary  (i.e., 
earth-fixed)  wave  of  s=2;  the  third  to  an  eastward  propagating  diurnal  tide  with  s=l.  We 
eliminate  from  consideration  (m,T)  since  this  represents  a  migrating  diurnal  tide 
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Figure  24.  Superposition  of  spectra  in  different  latitude  bands, 
corresponding  to  the  density  residuals  depicted  in  Figure  17  for  July,  1983. 


which  cannot  be  detected  from  a  single  fixed-local  time  measurement  (in  other  words,  it 
does  not  meet  the  Nyquist  criterion).  The  only  plausible  choice  appears  to  be  the  quasi¬ 
stationary  wave.  For  ks  =  4,  we  have  (m,T)  =  (4,  oo),  (0, 1/4),  (5,-1),  etc.  Again,  the 
migrating  semidiurnal  tide  (2, 1  /2)  is  eliminated  from  consideration,  and  there  are  no 
known  reasons  to  expect  an  eastward-propagating  wave  or  a  zonally-symmetric  6- hour  tide. 
Our  tentative  conclusion  is  that  the  structures  seen  in  Figure  23,  at  least  those  consistent 
with  apparent  zonal  wavenumbers  of  4,  2,  3,  5,  and  1  in  Figure  24,  correspond  to  earth-fixed 
features.  It  is  noteworthy  that  the  most  prominent  wavenumber,  s  =  4,  also  corresponds 
to  the  most  prominent  wavenumber  in  the  decomposition  of  land-sea  surfaces.  This  may 
suggest  some  relationship  with  lower-atmosphere  physical  processes.  The  thermospheric 
density  amplitude  variations  associated  with  these  features  are  about  ±2-4%. 

We  are  led  to  the  following  conclusion:  Longitudinal  wave  structures  of 
order  ±2-4%  exist  in  thermospheric  density,  with  apparent  earth-fixed  zonal 
wavenumbers  ranging  between  1  and  5.  Zonal  wavenumber  4  appears  to  play  an 
important  role.  Analyses  of  additional  geophysical  periods,  perhaps  using  more 
sophisticated  spectral  analysis  techniques,  may  lead  to  improved  quantification 
and  interpretation  of  this  phenomenon. 
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